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We calculate angle-dependent heat capacity in a low magnetic field range on the basis of Kramer- 
Pesch approximation together with an electronic structure obtained by first-principles calculations 
to determine a superconducting gap function of (TMTSF)2C104 through its comparisons with ex- 
periments. The present comparative studies reveal that a nodal d-wave gap function consistently 
explains the experimental results for (TMTSF)2C104. Especially, it is emphasized that the observed 
unusual axis-asymmetry of the angle-dependence eliminates the possibility of s-wave and node-less 
d-wave functions. It is also found that the directional ordering of CIO4 anions does not have any 
significant effects on the Fermi surface structure contrary to the previous modelings since the two 
Fermi surfaces obtained by the band calculations almost cross within the present full accuracy in 
first-principles calculations. 

PACS numbers: 71.15.Mb, 74.25. Jb, 74.20.Rp, 74.25. Bt 



I. INTRODUCTION 

The organic superconductors, (TMTSF)2X's have at- 
tracted considerable attention because of their rich va- 
riety of phases including superconducting states in spite 
of rather simple quasi-one-dimensional (QID) electronic 
structures. For example, (TMTSF)2PF6 exhibits a pres- 
sure induced superconductivity whose pairing symme- 
try has been expected to be spin-triplet by Knight-shift 
measurementsi. Meanwhile, (TMTSF)2C104 IS an am- 
bient pressure superconductor, whose Knight-shift mea- 
surements suggest a conventional spin-singlet pairing^. 
Since the singlet superconducting phase is in close prox- 
imity to spin density wave (SDW) states, spin-fluctuation 
has been proposed as a candidate of the pairing glue. 
However, the pairing mechanism is still far from an es- 
tablished settlement for (TMTSF)2J^ in spite of several 
theoretical studies^"—. 

Synchronous to the theoretical struggle, a num- 
ber of experimental studies have been also made 
to clarify the superconducting pairing symmetry for 
(TMTSF)2C104. The superconducting transition tem- 
perature for (TMTSF)2(C104)i-x(Re04):c is suppressed 
by a tiny amount of non-magnetic impurities. Such a 
result may be evidence of the presence of nodes in the su- 
perconducting gap function. Moreover, Takigawa et al. 
found that the nuclear magnetic relaxation rate lacks a 
coherence peak below Tc together with a low temperature 
power-law behavior (1/Ti ~ T^)^. These results suggest 
that the gap function has line-nodes at which the sign of 
the gap function changes. On the other hand, Belin and 
Behnia showed that the thermal conductivity rapidly de- 
creases below Tc- The data leads to the absence of any 
nodal structuresi^. 

Recently, a technique rotating the applied magnetic 
field has been incorporated in thermal measurements 



such as heat capacity and thermal conductivity measure- 
ments to probe the gap structure including positions of 
gap nodes in more detailsi^ii^. The angle-resolved heat 
capacity measurement is one of such advanced measure- 
ments, in which one detects details of the gap structure, 
especially locations of its nodes by measuring the oscil- 
lation of the heat capacity C{H) with respect to the 
applied magnetic field H direction. The rotational de- 
pendence on the thermal conductivity is also a power- 
ful tool to examine gap structures similarly. In this pa- 
per, we propose Kramer-Pesch approximation (KPA) to- 
gether with first-principles electronic-structure calcula- 
tions as a new theoretical tool to analyze the advanced 
angle-dependent-measurement datapi^. KPA significantly 
exceeds the accuracy level of the previous Doppler shift 
approximation. We have calculated the density of states 
(DOS's) around a vortex by using anisotropic Fermi 
surfacesii"— obtained by first-principles calculations and 
compare calculated angle dependences of heat capacity 
with the experimental data. 

Very recently, Yonezawa et al. reported that the 
oscillation curve of the heat capacity C{H)/T for 
(TMTSF)2C104 obtained by rotating the field becomes 
asymmetric with respect to the crystalline a-axis in the 
low temperature and low field range. In addition, kink 
structures are observed near this a-axis direction. They 
claimed, based on the Doppler shift analysis, that the 
asymmetries and the kink structures are clear evidence 
for the presence of line-nodes. 

On the other hand, a structural phase transition cor- 
responding to an anion ordering was observed at 24 K 
for (TMTSF)2C10^-^. Shimahara proposed a nodeless 
fully-gapped rf-wave superconductivity associated with 
the anion order^. Afterwards, several theoretical studies 
have suggested effective models with or without the anion 
ordering to discuss the pairing mechanism. '^^''■^ However, 
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there has been no first-principles calculation taking ac- 
count of the anion ordering. Here, we emphasize that 
a trustful band calculation resolves effects of the anion 
ordering on the electronic structure. 

The present paper has two objectives. The first 
one is to analyze the angle-resolved heat capacity for 
(TMTSF)2C104 with Fermi surfaces obtained by first- 
principles calculations. We identify the superconduct- 
ing gap function by using an expression for the heat ca- 
pacity on the basis of KPA in a low temperature and 
low magnetic field range. The KPA-based expression is 
applicable to various other unconventional superconduc- 
tors. We calculate the heat capacity by assuming the 
s-wave, nodeless d-wave and nodal d-wave gap functions 
and compare the results with the measurement data. The 
second one is to examine the effects of the anion order- 
ing on the electronic structure of (TMTSF)2C104. We 
clarify the electronic structure, especially Fermi surfaces 
through first-principles calculations using the measured 
structure parameters in the anion-ordered state. 

The rest of this paper is organized as follows. The 
quasiclassical approximation to describe the supercon- 
ducting state is briefly introduced in Sec. II. Then, the 
KPA calculation scheme based on the quasiclassical ap- 
proach is presented in Sec. III. We derive a vortex so- 
lution using KPA and present an expression of the heat 
capacity around a vortex. The electronic structure of 
(TMTSF)2C104 by first-principles calculations is given in 
Sec. IV, in which we display the band structure and Fermi 
surfaces. The calculation results on the angle-dependent 
heat capacity are shown in Sec. V. The discussion and 
conclusion are, respectively, given in Sec. VI. and VII. 



II. QUASICLASSICAL THEORY OF 
SUPERCONDUCTIVITY 

In many BCS superconductors, the gap-amplitude 
is much smaller than the Fermi energy, |A| <C Ep. 
In this case, one can properly use a quasiclassical 
approximation- ^-^^ . We consider the quasiclassical 
Green's function g that has the matrix elements in the 
Nambu (particle-hole) space as 

giz,r,kp)^(^f^ /J, (1) 

which is a 2 X 2 matrix in the Nambu space and is a 
function of complex frequency z , Fermi wave- vector k p , 
and point r in real space. We set /i = fee = 1 through 
this paper. The equation of motion for g is written as 

-it;F(feF) ■ Vg= [zf3 - A(r,fcF),ff] , (2) 

with Fermi velocity k-p and the commutator [a, h] = ab- 
ba supplemented by the normalization condition 

f = -vr^i. (3) 



Here, A is given by 

We neglect the vector potential by confining ourselves 
in type II limit. Setting z = e -\- irj with infinitesimal 
positive T], we obtain the retarded quasiclassical Green's 
function g^. In this paper, we use a special parameteri- 
zation form of the quasiclassical Green's function to solve 
Eq. ^i^-^ The solution g of Eq. can be written as 

_ _ -ITT / 1 - a6 2ia \ , ^ 

where a and b are the solutions of the following Riccati 
differential equations: 

vp{kp) ■ V a = 2iza - a^ A* + (6) 
i;F(fcF)- Vfe=-2iz6 + fe2A-A*. (7) 

In the parameterization Eq. ([5]), the normalization con- 
dition Eq. ([3]) is automatically satisfied. '^^ 

Since Eqs. ^ and ([7]) contain V only through vp{kp) ■ 
V, they are reduced to a one-dimensional problem on a 
straight line, the direction of which is given by that of 
the Fermi velocity vp{kp). We consider a single vortex 
along the zm axis. Because of the translational symmetry 
along the axis, the pair potential A does not depend 
on zm in the Riccati equations ^ and ([7]) , and hence the 
Riccati equations can be rewritten as 

d 

vp±(kp)—a — 2iza — a^A*(s, y, kp) + A(s, y, kp), (8) 
OS 

d 

vp±{kp)—b^-2izb + b^A{s,y,kp)-A*{s,y,kp), 
OS 

(9) 

where vp±{kp) is the amplitude of the vector vp±{kp) 
perpendicular to the zm axis by projecting the Fermi ve- 
locity v{kp) and the coordinate s (y) is along the direc- 
tion parallel (perpendicular) to vp±{kp). For simplicity, 
we solve the Riccati differential equations under a given 
form of the pair function. The density of states is given 
by 

7V(6) = (Kr,e))sp, (10) 
-(r,6)---/^^Im(ff^). (11) 

Here, (••■)sp = Jo° i^dr J^^ da/{Trrl) is the real-space 
average around a vortex where Ta/^o — \J Hc2 / H [Hc2 = 
^o/iT^^o), $0 — T^fa-iH] and dSp is the Fermi-surface area 
element. By using N{e), the low-temperature specific 
heat is given as 

Cm^r de 6^ N{e) 
T T T2 cosh^ {4^y ^ ' 



3 



III. KRAMER-PESCH APPROXIMATION 

We introduce Kramer-Pesch approximation (KPA) as 
an efficient method to analyze the angle-resolved exper- 
iments. We have shown that KPA gives the zero-energy 
density of states around a vortex consistent quantita- 
tively with results of direct numerical calculations^^. In 
addition, the computational time required for KPA is 
almost the same as that for the Doppler Shift method, 
which is significantly less than that in direct numerical 
calculations. Furthermore, it is emphasized that KPA 
can calculate the density of states even in complicated 
Fermi surfaces without any heavy numerical computa- 
tions. So far, we have actually examined various uncon- 
ventional superconductors with the use of KPAi^"— i^i^^ 

In works based on KPA, there has been a different way 
in theoretical treatments on the vortex core. Mel'nikov 
et al. presented an analytical solution describing the 
anomalous branches in a single vortex with arbitrary 
winding numbers by generalizing the Caroli-de Gennes- 
Matricon approach'^'^. They also demonstrated that the 
analytical solution on a single vortex is valid even in a 
higher energy range near the gap-amplitude. Therefore, 
we incorporate the Mel'nikov's method to calculate the 
heat capacity. The Mel'nikov's method can be regarded 
as a perturbation with respect to both energy and imag- 
inary part of the pair-function in the Riccati formalism. 

Now, let us show the present scheme. First, we briefly 
mention the Doppler Shift method in the Riccati formal- 
ism for comparison. We separate the pair-potential A 
into the amplitude and the phase $(s,j/, fcp) as 

A(s,y,feF) = |A(s,2;,fcF)|e^*(^'^''=-^ (13) 

Introducing a = exp(i$)a and b = exp(— i<I>)6, the Ric- 
cati equation is written as 

d ( d \ 

VF±{kF)-g^a^il2z~VFA^{kF)-^^ja~a^\A\ + \A\, 



(14) 



VF±ikF)-^b : 

as 



-I I 20-WFi(fcF)^$ ) 6 + fo'|A| - lAI 



(15) 



Assuming da/ds = db/ds = 0, the equations can be ex- 
actly solved in an analytical way. The solution is equiv- 
alent to that in the bulk region by replacing the energy 
z with the Doppler shifted energy z — {vF±/2){d^/ds). 
The Doppler Shift method is an approximation neglect- 
ing the spatial variation of \a\ = \a\. Then, it breaks 
down near a vortex core^^i^i. 

Next, we derive the vortex solution by using KPA. We 
write down a pair-potential around a vortex in the fol- 
lowing form, 

A(s,y,fcF) = f{s,y)A^d{kF)e'''-, (16) 



= /(.,2;)Ao.d(fcF)^^e^^"('=^) 



(17) 



where f{s,y) describes the spatial variation of the pair- 
potential. Then, /(O) = 0, limj—j-oo /(?") = 1, and Aoo 
is a pair-potential in the bulk region. 9r denotes an an- 
gle around a vortex and 9y does a direction of the pro- 
jected Fermi velocity VF±{kF)^'^. Introducing the vari- 
ables written as 



a = ae " , 
b = 6e-'^' 
A = Ae'^" 



the Riccati equations are rewritten as 

VF±ikF)^a = 2iza - a^K* + A, 
ds 



VF^{kF)^b 

ds 



'2iz6 + 6 A 



(18) 
(19) 
(20) 

(21) 
(22) 



In KPA with the Riccati formalism as the previous 
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we expand a and b in Eqs. (j2ip and 



paper- 

with respect to the impact parameter y and the complex 
frequency z. In this paper, we expand these variables 
with respect to the imaginary part of the pair-function 
A, in stead of y, and the complex frequency z on the 
basis of the Mel'nikov's method. Then, Ar and Aj are 
defined as 



AR = ReA = /(s,y)Aood(fcF) 
Ai =ImA /(s,y)Aood(feF) 



(23) 
(24) 



Following Refs. |25| and |32|, we eventually obtain a and b 
as 



a = ao + ai + 0(z^, A^, zAi), 
b = bo + bi + Oiz^,AlzAi), 



with 



ao = -sign (rf(fcF)), 
bo = sign {d{kF)), 



(25) 
(26) 

(27) 
(28) 



ai(s) 



his) 



VF±ikF) 
VFA_{kF) 



(2moz + 2iAi(s'))e""^"'^ds' (29) 

OO 

(-2i6oz + 2iAi(s'))e~"^"'^ds'. 



(30) 



Together with the help of the following function 



u{s) 



A^f{s',y)^==ds', (31) 

V '5 + y 



VF±ikF 

the quasiclassical Green's function is then written as 
-27ri 
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aibo + aob_ 

7rWF±(feF) 



■M, 



C(2/,feF) z-E{y,kF) 



M, 



(32) 
(33) 
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with 

M 

C{y,kF) 



1 iflo 
-ibo -1 



_ M(fcF)|Ac 



(34) 
(35) 

(36) 



The quasiclassical Green's function has the pole at z = 
E{ii,k-p), which is regarded as energy of quasi-particle. 
Substituting Eq. ^ into Eqs. ^ and ^ and setting 
z = e + «77, the density of states is given as 



N{e) 



dSv UF±(feF) _„ 



27r2|t;F|C(y,fcF: 



e-«(^)^(e_iJ(y,fcF)))^p. 

(37) 



Thus, we obtain the heat capacity in the clean-hmit (ry 
0) with the use of the KPA written as 



C(T) 
T 



4 s) 



2^2|^p|T3 Civ.kv) ^Q3^^2 (^ E{yM^) \^ /sp' 

(38) 

On the other hand, we set the spatial variation of the 
pair-potential /(s, y) as 



(39) 



With the use of this function, one can integrate Eqs. ([35 
and dMl): 



C(y, kp) = 2^2/2 + QKiiroiy, fcp)), 



E{y,kF) = |d(fcF)|A 
with 

ro{y,kF 



y Ko{ro{y,kF)) 



y7T^Xi(ro(2/,fcF))^ 

2|d(fcF)|A^ r— 



(40) 
(41) 

(42) 



where, the function K„(x) is the modified Bessel function 
of the second kind. 




FIG. 1. (Color online) Crystal structure for (TMTSF)2C104. 



as accuracy. The lattice constants and atomic inner- 
coordinates refer to a measurement report, Ref. Issl . The 
calculation self-consistent loops to obtain a converged 
electronic structure are repeated until the total energy 
difference becomes smaller than 10 eV. In the loops, 
/c-points are taken as 10 x 5 x 5, and the energy cut-off is 
set to be 500 eV. Once the electron density is obtained 
after the convergence, the energy bands are again calcu- 
lated on finer fc-points as 49 x 23 x 25 in order to determine 
the Fermi surfaces and the Fermi velocities accurately as 
much as possible. 

Since the employed structural parameters are mea- 
sured at 7K, the data reflects the orientational ordering 
of the tetrahedral CIO4 anions in the crystal structure of 
(TMTSF)4C104. The ordered structure is displayed in 
Fig. [TJ We obtain the band structure and Fermi surface 
for the structural parameters as shown in Fig. [5] and [31 
From Fig. [31 it is found that the two Fermi surfaces al- 
most cross each other since the anion ordering gap is too 
small to resolve it in the standard scale. The tiny gap can 
be distinguished only by an enlarged scope as the inset 
of Fig. 3. This result clearly suggests that the direction 
of the CIO4 anion ordering does not have any significant 
effect on the Fermi surface structure contrary to the pre- 
vious theoretical expectations. Thus, we would like to 
point out that any theoretical modelings originated from 
the anion ordering are unlikely to consistently explain 
recent advanced experimental data. 



RESULTS 



IV. BAND STRUCTURE 

Now, let us display the electronic structure. In order to 
calculate the band structure, we employ a first-principles 
density-functional-calculation package VASP— . Among 
available options for the band calculations, we adopt 
GGA exchange-correlation energy^^ and PAW method'^^ 
due to their excellent computational performance as well 



We study the angle-resolved heat capacity, in which 
the applied magnetic field is rotated inside the basal a- 
6-plane. Three pairing symmetries are employed to test 
their matching with the angle-resolved experimental re- 
sults. For simplicity, we assume that the vortex core is 
cylindrically isotropic and anisotropy of the critical mag- 
netic field is not present, i.e., Hc2{4>) ~ Hc2- On the anal- 
ysis of the vortex core excitation, we set a spatial cutoff 
length = 5^0 J which is comparable to the neighboring 
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s-wave gap function 



FIG. 2. (Color online) Band structure for (TMTSF)2C104. 




FIG. 3. (Color online) Fermi surfaces for (TMTSF)2C104. 
Inset: Closeup of the Fermi-surface crossing, sliced at = 
0. The mesh describes fc-points actually used in the Fermi- 
surface calculation. 



vortex distance as the magnetic field H ^ Hc2/25. We 
take the x-axis (y-axis) parallel (perpendicular) to the 
a-axis. It is also noted thaty-axis is parallel to the b'- 
axis introduced by the Ref. |39| on a-&-plane and z-axis is 
perpendicular to the a-6-plane. 



A. s-wave gap function 

First, we examine a possibility of an isotropic s-wave 
gap function. In this case, the oscillation pattern of the 
angle resolved heat capacity suffers only the Fermi sur- 
face anisotropy. As shown in Fig. 21 the heat capac- 
ity curve monotonically oscillates with the angle and 
shows the minima at = reflecting the Fermi surface 
anisotropy. The minima appear as the magnetic field di- 
rection is parallel to the a-axis. This oscillatory pattern 
is as one expects, but inconsistent with the latest mea- 
surement data in details. In terms of the minima, we note 
that the Doppler shift method can not resolve even these 
minima, since the Doppler shift method can not describe 
the Fermi surface anisotropy in fully-gapped supercon- 
ductors. Meanwhile, in the present scheme using KPA, 
the momentum fcp dependent kernel of the heat capacity 
(the integrand in Eq. ([55)) ) vanishes as the magnetic field 
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FIG. 4. (Color online) Angular dependence of the heat ca- 
pacity rotating magnetic fields on a-b' plane in the case of the 
s-wave gap function. 



is directed parallel to the Fermi velocity v{kp), since the 
projected Fermi velocity wfj^C^f) then becomes zero. In 
the case of (TMTSF)2C104, the Fermi velocity is almost 
parallel to the a-axis on the whole Fermi surfaces because 
of the quasi-one-dimensionality. 



B. Nodeless d-wave gap function 

Next, we check a possibility of the nodeless rf-wave 
gap function proposed by Shimaharaj^ We employ the 
nodeless d-wave gap function expressed as (see Fig. [5]) 

A°(fc,, ky, fc,) = ao/°(fc„ fey + bo), (43) 
A\k,, ky, k,) = -{aif\k,,ky) + 6i), (44) 

with 

ky) = cos (^ky - ^o)^^ - ^ cos (^{ky - ^o)^^ 

(46) 
(47) 



f^{kx,ky) = cos ^(fcy - ko)^-^ 
kf) = —0.031 sign {k^). 



where, ao,i and &o,i denote normalization factors (ao = 
10/26, bo '= -17/26, ai = 10/23 and bi = -3/23), gy is 
the second element of the reciprocal lattice vector defined 
by gy — 27t/ (6sin7), b denotes the crystal axis b — 15.356 
A , 7 = 68.92°, and /cq denotes the intersection of two 
Fermi surfaces. As shown in Fig. [6l a slightly con- 
cave but almost fiat-like curve showing the minimum at 
the center lies around the angles —20° < cj) < 20°. These 
symmetric kink-like curvature including the unclear min- 
imum is due to anisotropy of the Fermi surface and the 
gap function. Close to ky ~ kQ, the Fermi velocity con- 
tinuously changes its direction around the a-axis, and the 
amplitude of the gap function gives the minimum with 
no change of the gap sign as shown in Fig. [S] Therefore, 
the angle variation of the magnetic field around a-axis 
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FIG. 5. (Color online) Gap functions on the Fermi surfaces 
in the case of the nodeless d-wave gap function. 



nodeless d-wave gap function 
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FIG. 6. (Color online) Angular dependence of the heat ca- 
pacity rotating magnetic fields on a-b' plane in the case of the 
nodeless d-wave gap function. 




FIG. 7. (Color online) Schematic figure about the nodal line 
in the case of the nodal d-wave gap function I. 



The case I, whose nodes are on both the Fermi surfaces, 
is displayed in Fig. [71 while the case II (III) , whose nodes 
are only on the inner (outer) Fermi surface. 

As shown in Fig. [HI the asymmetric behavior with re- 
spect to the a-axis {(j> = 0°) direction is found in the 
case of the nodal gap function II. Paying attention on 
the curve near = 0° as shown in Fig. [9l the kink like 
structures are observed at around (j) — ±15° in an anti- 
symmetric manner. This calculation using the nodal gap 
function II is mostly consistent with the experimental 
resultsjii^ These results suggest that the superconducting 
gap nodes He on the inner Fermi surface in the 0-th band. 



VI. DISCUSSION 



(-20° <(j)< 20°) is alnrost lost, i.e., DOS's of the quasi- 
particles with the small gap do not almost change with 
the angle. These results are inconsistent with the mea- 
surement data. 



C. Nodal d-wave gap function 

Finally, we examine nodal c?-wave gap functions. The 
trial nodal d-wave gap functions are classified into three 
types expressed as 

{f{kx,ky) (easel) 
f{k,,ky) (case II) , (48) 
1 (case III) 

{f{kx,ky) (easel) 
1 (case II) , (49) 

f{kx,ky) (casein) 

with 

f(kx, ky) = cos {^y-^ ■ (50) 



First, we mention the electronic structure including the 
Fermi surfaces. In order to compare the present one with 
the previous tight-binding models^iiiS, we construct the 
tight-binding model based on our band calculation tak- 
ing account of the anion ordering. As shown in Fig. IIOI 
the obtained tight-binding model is equivalent with that 
by the band calculation. Following the Ref. [2l| in terms 
of the notation for the transfer integrals, the values we 
obtain are is,^ = 280 meV, is,^ = 247 meV, ts^^ = 269 
meV, ts2B = 248 meV, t/^ = -47.0 meV, ti^ = -57.9 
meV, t/3 = 48.0 meV, t/, = -10.2 meV, t/^, = 63.3 meV, 
and t/g = 3.98 meV. These values are much more close to 
the values obtained by the band calculation without the 
anion ordering^ than those obtained by phenomenolog- 
ically taking account of the anion ordering^?- Therefore, 
the Fermi surfaces shown in Fig. [3|are found to be qual- 
itatively different from those obtained by the previous 
tight-binding model considering the anion ordering. 

Next, we discuss the assumption used in the gap ex- 
amination, i.e., the excitation structure around vortex 
core is isotropic. In the QID superconductors, the su- 
perconducting gap amplitude may spatially vary around 
a vortex. However, low-lying quasiparticle excitations 
around a vortex are usually less affected by the variation 
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nodal d-wave gap function I 
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nodal d-wave gap function II 
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nodal d-wave gap function III 




-45 45 90 135 180 
(|) [degree] 



FIG. 8. (Color online) Angular dependence of the heat ca- 
pacity rotating magnetic fields on a-b' plane in the case of the 
nodal d-wave gap function I, II, and III. 
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FIG. 10. (Color online) tight-binding fitted band structure 
(solid curves), the squares denote the first- principle band 
structure shown in Fig. [2] 



of the gap amplitude than that of the superconducting 
phasei^2iiL Our expression on the heat capacity shown in 
Eq. (p8|) can contain an anisotropic vortex core through 
/(s, y). Therefore, it is a next-step checkpoint for the ex- 
perimental consistency to investigate effects of the vortex 
core including the anisotropy. 

Finally, we discuss the origin of the asymmetric be- 
havior in the curves of the angle-resolved heat capacity. 
We calculate the partial heat capacity from each Fermi 
surface in the case of the nodal d-wave gap function I. 
As shown in Fig. [TTl the kink structures of the angle- 
dependent partial heat capacity are different on each 
Fermi surface at T = 0.05A. This difference originates 
from that in the Fermi velocity of the nodal quasiparti- 
cles on each Fermi surface. One finds that the asymmet- 
ric behaviors of the partial heat capacity are observed 
only in the case with the gap-nodes, since the distribu- 
tion of the direction of the Fermi velocity does not have 
the strong asymmetry on the whole Fermi surfaces of 
(TMTSF)2C104. Therefore, the asymmetric kink struc- 
ture of the curves measured in the angle-resolved heat 
capacity is a clear evidence that the gap function has 
nodes in (TMTSF)2C104. 



nodal d-wave gap function II 
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FIG. 9. (Color online) Angular dependence of the heat ca- 
pacity rotating magnetic fields on a-b' plane in the case of the 
nodal d-wave gap function II. 



VII. CONCLUSION 

In order to examine the effects of the anion order- 
ing and resolve the superconducting gap function in 
the organic superconductor (TMTSF)2C104, we per- 
formed first-principles calculations and developed the 
quasi-classical theory, respectively. The first-principles 
calculation revealed that the anion ordering does not 
have any important role on the Fermi surface shapes in 
contrast that the gap opening around the crossing point 
was previously expected as a consequence of the order- 
ing. The present calculation partly excludes the previ- 
ous modeling based on the intuitive expectation. On the 
other hand, using Kramer-Pesch approximation on the 
single vortex core excitation together with the Fermi sur- 
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nodal d-wave gap function I at T = 0.05 A 
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FIG. 11. (Color online) Angular dependence of the partial 
heat capacity rotating magnetic fields on a-h' plane in the 
case of the nodal d-wave gap function 1 at T = 0.05A. Ni is 
the partial heat capacity on the i-th band. 



faces obtained by the first-principles calculations, we con- 
structed the formula calculating the angle-resolved heat 
capacity in the low-field range and compare the angle 
dependence obtained from various gap function models 
with the experimental results. Consequently, we showed 
that the nodal d-wave gap function consistently explains 
the experimental results. Especially, it should be em- 
phasized that only the nodal d-wave gap function can 
reproduce the axis asymmetry of the angle dependence. 
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